-- Create a frequency vs. time  spectrogram image.
-- Requires Xfun "fft"

-- example parabolic freq sweep
sweep(t) = cos(2**Radians*(.1+.015*t^2)*t)

ndata=512
dt=1/30
signal[i] = sweep(dt*(i-1)) dim[ndata]  -- sample it

Xmin=0;  Xmax=dt*ndata;  Xlabel="sec"
plotline signal

-- use a sliding data window to see power spectrum changes
-- "window" is the number of data points in the window
-- "slide" is the number of data points to move the window for each spectra

spect(signal,window,slide)[i,j] = i0:=(i-1)*slide,
      pwrspec(signal[i0+1:i0+window]) when j=1,,
      transform[j]    dim[count(signal)/slide,window/2+1]

window=64    -- must be a power of 2
slide:=trunc(window/4):
-- show window and slide time on signal plot
plotline {{0,window*dt},{-.93,-.93}} --{{x1,x2},{y1,y2}}
plotline {{slide,slide+window}*dt,{-.9,-.9}}

newaxis
Zmin=0; Zmax=.3
Ymin:=-1/(window*dt):; Ymax:=1/(2*dt):; Ylabel:="  Hz":
image  spect(signal,window,slide)

-- Note that when the window is moved by less than its width the last few spectra will be undefined because there is not a full window of data available.
